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We review recent progress of the HAL QCD method which was recently proposed to in- 
vestigate hadron interactions in lattice QCD. The strategy to extract the energy-independent 
non-local potential in lattice QCD is explained in detail. The method is applied to study 
nucleon-nucleon, nucleon-hyperon, hyperon-hyperon and meson-baryon interactions. Several 
(-H . extensions of the method are also discussed. 

§1. Introduction 

00 ■ Oire of the ultimate goals in nuclear physics is to describe hadronic many-body 

problems on the basis of the hadronic S-matrices calculated from first principle QCD. 
■ In particular, the nuclear forces are the most fundamental quantities: Once they are 

obtained from QCD, one can solve finite nuclei, hypernuclei, nuclear matter and 
hyperon matter by employing various many-body techniques developed in nuclear 
CN ■ physics. 

Phenomenological nucleon-nucleon (NN) potentials, which are designed to re- 
produce a large number of proton-proton and neutron-proton scattering data as well 
as deuteron properties have been constructed in 90's and are called high-precision 
AW potentials. Some of the examples are shown in Fig. 1, which reflect character- 
istic features of the AW interaction for different values of the relative distance r as 
reviewed in l)-5): 

The long range part of the AW force (r > 2 fm) is dominated by one-pion exchange 
originally introduced by Yukawa. 10 ) Since the pion is the Nambu-Goldstone boson 
associated with the spontaneous breaking of chiral symmetry, it couples to the nu- 
cleon's spin-isospin density and leads to not only the central force but also the tensor 
force. The medium range part (1 fm < r < 2 fm) of the AW force receives signifi- 
cant contributions from two-pion (irir) exchange 11 ) and/or heavy meson (p, oj, and 
a) exchanges. In particular, the spin-isospin independent attraction of about 50 - 
100 MeV in this region plays an essential role to bind the atomic nuclei and nuclear 
matter. The short range part (r < 1 fm) of the AW force is best described by a 
phenomenological repulsive core introduced by Jastrow. 12 ) 
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Fig. 1. Three examples of the modern NN potential in 1 So (spin-singlet and S- wave) channel: 
Bonn, 6 ' Reid93 7) and Argonne v 18 . 8) Taken from Ref. 9). 



The nuclear saturation, the nuclear shell structure, the nuclear superfluidity 
and the structure of neutron stars are all related to the properties of the nuclear 
force. 13 )~ 15 ) Furthermore, the hyperon-nucleon (YN) and hyperon-hyperon (YY) 
forces, whose information is still quite limited experimentally, are crucial to under- 
stand the structure of hypernuclei and the core of the neutron stars. The three- 
nucleon forces (and the three-baryon forces in general) are also important to under- 
stand the binding energies of finite nuclei and the equation of state of dense hadronic 
matter. 

It has been a long-standing challenge in theoretical particle and nuclear physics 
to extract the hadron-hadron interactions from first principle. A framework suitable 
for such a purpose in lattice QCD was first proposed by Liischer: 16 ) For two hadrons 
in a finite box with the size L x L x L under periodic boundary conditions, an 
exact relation between the energy spectra in the box and the elastic scattering phase 
shift at these energies has been derived. If the range of the hadron interaction R is 
sufficiently smaller than the size of the box R < L/2, the behavior of the two-particle 
Nambu-Bethe-Salpeter (NBS) wave function (p(r) in the interval R < |r| < L/2 is 
sufficient to relate the phase shift and the two-particle spectrum. This Liischer's 
finite volume method bypasses the difficulty to treat the real-time scattering process 
on the Euclidean lattice. Furthermore, it utilizes the finiteness of the lattice box 
effectively to extract the information of the on-shell scattering matrix and the phase 
shift. 

A closely related but a new approach to the hadron interactions from lattice QCD 
has been proposed recently by three of the present authors 9 )' 17 )' 18 ) and has been 
developed extensively by the HAL QCD Collaboration. (Therefore the approach 
is now called the HAL QCD method.) Its starting point is the same NBS wave 
function (p(r) as discussed in Ref. 16). Instead of looking at the wave function 
outside the range of the interaction, the authors consider the internal region \r\ < R 
and define an integral kernel (or the non-local "potential" in short) U(r,r') from 
tp(r) so that it obeys the Schrodinger type equation in a finite box. This potential 
can be shown to be energy-independent by construction. Since U(r,r') for strong 
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interactions is localized in its spatial coordinates due to confinement of quarks and 
gluons, it receives only weak finite volume effect in a large box. Therefore, once U 
is determined and is appropriately extrapolated to L — > oo, one may simply use the 
Schrodinger type equation in infinite space to calculate the scattering phase shifts 
and bound state spectra to compare the results with experimental data. Since U is 
a smooth function of the quark masses, it is relatively easy to handle on the lattice. 
This is in sharp contrast to the scattering length, which shows a singular behavior 
in the quark mass corresponding to the formation of the hadronic bound state. A 
further advantage of the HAL QCD method is that it can be generalized directly to 
the many-body forces and also to the case of inelastic scattering. 

Studying structure of S (strangeness) = — 1 and S = —2 hypernuclei is one of the 
key challenges in modern nuclear physics. Also, the central core of the neutron stars 
will have hyperonic matter if the neutron beta-decays to hyperons become possible 
at high density. The hyperon-nucleon (YN) and hyperon-hyperon (YY) interactions 
are crucial to determine the level structures of hypernuclei as well as onset-density 
of hyperonic matter in neutron stars. 19 - > By generalizing the NN scattering in the 
flavor SU(2) space to the baryon-baryon (BB) scatterings in the flavor SU(3) space, 
the HAL QCD method can give the YN and YY potentials as natural extension 
of the NN potentials. Such extension is also useful for identifying the origin of the 
short-range repulsive core of the NN potential and for studying possible S = —2 
six-quark state such as the iJ-dibaryon. 

In this article, we review the basic ideas and recent progress of the HAL QCD 
method to hadron interactions. (As for the Liischer's finite volume method, see a 
recent review Ref. 20).) In Sec. 2, the basic strategy to define the NN potential 
in QCD is explained. In Sec. 3, we introduce lattice formulations of the time- 
independent HAL QCD method originally proposed in Refs.9), 17), 18) as well as 
its time-dependent generalization. In Sec. 4, some recent results of lattice QCD cal- 
culations for the NN potential are given in both quenched and full QCD. Magnitude 
of the non-locality in U is also discussed in the section. In Sec. 5, the method is 
applied to the hyperon-nucleon interactions such as NS and NA systems. In Sec. 6, 
interactions between octet baryons are investigated in the flavor SU(3) limit, where 
up, down and strange quark masses are all equal. In Sec. 7, a generalization of the 
HAL QCD method to the case of inelastic scattering is given. In Sec. 8, we show 
results of the three-nucleon potential, especially its short distant structure. In Sec. 9, 
an application to the kaon-nucleon scattering is considered. Sec. 10 is devoted to 
summary and concluding remarks. 

§2. Defining the potential in QCD 

2.1. Nambu-Bethe-Salpeter (NBS) wave function 

A key quantity to define the baryon-baryon(£?f?) "potential" in QCD is the 
equal-time Nambu-Bethe-Salpeter wave function, 



<p w (x)e- wt = (0\T{B(r + x, t)B(r, t)}\2B, W, Sl s 2 ), 



(2-1) 
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where \2B, W, S1S2) is a QCD eigenstate for two baryons with equal mass mg, helicity 
Si and S2, total energy W = 2yJ k 2 + m 2 B , the relative momentum k, and the total 
momentum p (we take p = in this paper). Generalization to the unequal mass can 
be formulated in a similar manner. In the case of two nucleons, the local interpolating 
operator B(x) is taken as 

Ba(x) = ( ) = e abc (u T a (x)C l5 d b (x)) q c , a (x), q{x) = ( jgj ) , (2-2) 

where x = (x, t), a, b, c are the color indices, and a is the spinor index. The charge 
conjugation matrix is given by C = 7274, and p, n are proton and neutron operators 
while u, d denote up and down quark operators. Here tp w implicitly has two pairs of 
spinor-flavor indices from B a (r + x,t)Bp(r,t) as well as two helicity indices si and 

S2- 

The most important property of the above NBS wave function is as follows. If the 
total energy W lies below the threshold of meson production (i.e. W < 2ms + vtlm 
with the meson mass tum), it satisfies the Helmholtz equation with k = \k\ at 

r = \r\ — > 00, 

[fc 2 + V 2 ] f w (r) ~0. (2-3) 

Furthermore, the asymptotic behavior of the radial part of the NBS wave function 
for given orbital angular momentum L and total spin S reads 18 )' 21 ) 

w { LS) K sin(kr-L./2 + 5 LS (k)) el5Ls{k) _ 

kr 

Here Sis(k) is nothing but the phase shift obtained from the baryon-baryon S-matrix 
in QCD below the inelastic threshold. It should be remarked here that only the upper 
components of the spinor indices for the NBS wave function (a = 1, 2 and (3 = 1,2) 
are enough to reproduce all BB scattering phase shifts 5is(k) with L = 0, 1, 2, 3, ■ ■ ■ 
and S = 0, 1 (See Appendix A of Ref.18) for the precise expression of Eq.(2-4) and 
its relation to the S-matrix in QCD.) 

2.2. Non-local potential from the NBS wave function 

From the NBS wave function, we define a non-local potential through the rela- 

tion 9),17),18) 

(E k - H ) <p^(x) = J U af3n5 (x, y)^(y)d 3 y, f £ fc = |], H = ,(2-5) 

where U(x, y) is expected to be short-ranged because of absence of massless particle 
exchanges between two baryons. As mentioned in the previous subsection, it is 
enough to consider the upper spinor indices of a, /3, 7, 8: Then 16 components of 
Uafs-^s can be determined from 4 components of tp^ for 4 different combinations 
of (si,S2). Since the NBS wave function ip w is multiplicatively renormalized, the 
potential U(x,y) is finite and does not depend on the particular renormalization 
scheme. Note that, while Lorentz covariance is lost by using the equal-time NBS 
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wave function and Eq. (2-5) is written as a Schrodinger type equation, no non- 
relativistic approximation is employed here to define U(x,y). 

The non-local potential U(x, y) has been shown to be energy- independent. 17 ^' 
To see this, let V\h be the space spanned by the wave function at W < W t h = 
2m b + rriM'- Vth = {VcHW — ^th} where c represents quantum numbers of the 
NBS wave function other than energy W. Then the projection operator to Vth is 
given by 

P W *(X,V) = E E ^W^UWi, W 2 )^(y)i 

W U W 2 <W th C!,C 2 

= E E p c Wth (^;^i/) (2-6) 

W<W th c 

where N~^ C2 (Wi,W2) is defined as the inverse of the Hermitian operator 

JVcLC^Wb) = J ^WV^W^, W 1)2 < W th , (2-7) 

which satisfies 

in the restricted indices that < Wth- (We here assume that iV(Wi,W2) does 
not have zero eigenvalues in this restricted space.) 
Using these, the non-local potential is defined by 

U w *{x,y)= E E*W) [E k -H } ip ^(x)N-] C2 (W 1 ,W 2 ) lf ^(yy 

W 1; 2<W th Cl,C 2 

= E E ~ p c Wth ( w ; v)- (2-9) 

Then, it is easy to observe that the above non-local potential satisfies Eq.(2-5) at 
W < W th : 

Wi<W th ci 

= 9(W th - W) [E k - H ] <pY(x). (2-10) 

This non-local potential U(x, y) is energy independent by construction. It is also 
easy to see that we can make the potential local but energy-dependent. Similar 
trade-off between non-locality and energy-dependence has been also discussed long 
time ago in Ref.22) in a different context. Note however that the non-local potential 
U(x,y) which satisfied Eq. (2-5) at W < W^h is n °t unique. For example, one 
may add a term such as f(x)[l — P Wth (x,y)] with arbitrary functions f(x) to the 
non-local potential U(x,y) without affecting Eq.(2-5) at W < Wth- 
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We remark that One may define a non-local potential different from Eq.(2-9) as 
U°°(x,y) = Y, J2[E k -H }P™(W;x,y), (2-11) 

W<oo c 

which satisfies Eq. (2-5) for all W. This potential, however, becomes long-ranged, 
due to the presence of inelastic contributions above Wth- An extension of the HAL 
QCD method, which keeps the short-range nature of the potential while inelastic 
channels open, will be discussed in Sec. 7. 

The most general form of the Schrodinger type equation for the NBS wave func- 
tion has energy-dependent and non-local potential as shown in Ref.16). However, one 
can always remove its energy-dependence as demonstrated in the above derivation. 

2.3. Velocity expansion of the non-local potential 

If one knows NBS wave functions ip w for all W < W^, the non-local potential 
U can be constructed according to Eq. (2-9). In lattice QCD simulations in a finite 
box, however, only a limited number of wave functions at low energies (ground state 
and possibly a few low-lying excited states) can be obtained. In such a situation, it 
is useful to expand the non-local potential in terms of the velocity (derivative) with 
local coefficient functions; 23 ) 

U(x,y) = V(x,V)5 3 (x-y). (2-12) 

In the lowest few orders we have 

V(r, V) = Vb(r) + F CT (r)o-i • <x 2 + V T {r)S l2 + VLs(r)L • S +0(V 2 ), (2-13) 

" * ' " v ' 

LO NLO 

where r = \r\, (T{ is the Pauli-matrix acting on the spin index of the i-th baryon, 
S = {a i + <T2)/2 is the total spin, L = r x p is the angular momentum, and 



'12 



, (r • <ri)(r • cr 2 ) 



<7l • <T 2 (2-14) 



is the tensor operator. Each coefficient function is further decomposed into its flavor 
components. In the case of nucleons (i.e. Nf = 2 ), we have 

Vx(r) = V%(r) + VZ(r)T 1 -T 2 , X = 0, a, T, LS, • • • , (2-15) 

where Tj is the Pauli-matrix acting on the flavor index of the i-th nucleon. The form 
of the velocity expansion (2-13) agrees with the form determined by symmetries. 24 -* 
At the leading order of the velocity expansion, the local potential is given by 

V LO (r) = V (r) + V^(r)<n • cr 2 + V T (r)S 12 , (2-16) 

which is obtained from the NBS wave function at one value of W. Since S12 = for 
the spin-singlet state, for example, one has 



(2-17) 
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2.4. Remarks on the "scheme" -dependence of the potential 

We emphasize that the potential itself is not a physical observable, and is there- 
fore not a unique quantity in quantum mechanics and in field theory. In fact, the 
baryon-baryon potential in QCD depends on the choice of the interpolating baryon 
operator to define the NBS wave function. Among others, the local baryon operator 
used in HAL QCD method is a most convenient choice, since the reduction formula 
for composite particles can be derived in a simplest way for this choice. 25 )~ 27 ) 

Nevertheless, one may adopt other interpolating operators (such as higher di- 
mensional operators and non-local operators): Particular choice of the baryon op- 
erator and associated potential may be considered as a "scheme" to describe phys- 
ical observables such as the scattering phase shift and the binding energies. The 
potential, although being "scheme"-dependent, is still useful to understand physi- 
cal phenomena as we know well in quantum mechanics. The repulsive core of the 
nucleon-nucleon potential in the coordinate space, which is known to be the best 
way to summarize the NN scattering phase shift at high energies, is one of such 
examples.*) 

Among different schemes, good convergence of the velocity expansion is an im- 
portant check of the choice of our present scheme. Such a check can be carried out 
by examining the W dependence of the lower order potentials. For example, if we 
have (f Wn for n = 1,2,- • • N, we can determine the N — 1 unknown local functions 
of the velocity expansion in N different ways. The variation among N different de- 
terminations gives an estimate of the size of the higher order terms. Furthermore 
one of these higher order terms can be determined from ip Wn for n = 1, 2, ■ ■ ■ N. The 
convergence of the velocity expansion will be investigated explicitly in Sec. 4. 

The analysis in this section shows that the use of Schrodinger type equation 
with non-local potential is justified to describe the BB scattering in QCD. The key 
quantity is the NBS wave function, whose asymptotic behavior encodes phases of the 
S-matrix for the BB scattering. If the velocity expansion of the non-local potential 
is reasonably good at low energies, one can use the LO and NLO potentials to 
investigate various nuclear many-body problems. 

§3. Lattice formulation 

We now discuss procedures to extract the NBS wave function from lattice QCD 
simulations. For this purpose, we consider the correlation function on the lattice 
defined by 

F(r, t - t ) = (0\T{B(x + r, t)B(x, t)}J(t o )\0) (3-1) 

where J{to) is a source operator which creates two-baryon states. Inserting a com- 
plete set and considering baryon number conservation, we have 

F(r, t-to) = (0\T{B(x + r, t)B(x, t)} \2B, W n , s u s 2 )(2B, W n , s u s 2 \J(t )\0) 

n,si,S2 

*-* Although in a different sense of the "scheme", analogous situation in quantum field theory is 
the running coupling constant. It is scheme-dependent quantity but is quite useful to understand 
the high energy processes such as the deep inelastic scattering data. 



8 



S. Aoki et al. (HAL QCD Collaboration), 



+ ■■■= Y, A n ^ S2 ^{r)e- w ^-^ + • • • , (3-2) 

n,si,S2 

where A niS1)S2 = (2B, W n , s±, S2\J^(to)\0) an d ellipses represent contributions from 
inelastic states such as NNir, NNNN, etc. At large time separation (t — to) — > oo, 
we obtain 

lim F(r,t- t ) = A v Wo (r)e- Wo(t - u>) + O( e -^o(t-*o)) (3.3) 

(t— io)^oo 

where Wo is the lowest energy of BB states. Since the source dependent term Aq 
is just a multiplicative constant to the NBS wave function ip w °(r), the potential 
defined from ip w °(r) is manifestly source-independent. For this extraction of the 
wave function to work, the ground state saturation for F in Eq. (3-3) must be 
satisfied by taking large t — to- In practice, however, F becomes very noisy at large 
t — to. In Sec. 3.4, we will discuss more on this point. 

3.1. Choice of source operators 

We choose the source operator J to fix quantum numbers of \2B, W, s±, 82)- 
Since lattice QCD simulations are usually performed on a hyper-cubic lattice, the 
cubic transformation group £0(3, Z) instead of £0(3, R) is considered as the sym- 
metry of 3-dimensional space. Therefore the quantum number is classified in terms of 
the irreducible representation of SO(3, Z), which is denoted by A±, A2, E, T\ and T2 
whose dimensions are 1, 1, 2, 3 and 3, respectively. Relation of irreducible representa- 
tions between S0(3, Z) and 50(3, R) is given in Table I for L < 6, where L denotes 
the angular momentum for the irreducible representation of 50(3, R). For example, 
the source operator J' (to) in the A\ representation with positive parity generates 
states with L = 0,4, 6, ••• at t = to, while the operator in the T\ representation 
with negative parity produces states with L = 1,3,5, ••• . For two octet-baryons, 
the total spin S becomes 1/2 ® 1/2 = 1 © 0, which corresponds to T\(S = 1) and 
Ai(S = 0) of S*0(3, Z). The total representation J for a two baryon system is thus 
determined by the product R\ (g> R2, where R\ = A\, A2, E, T\, T2 for the orbital 
"angular momentum" while R2 = A±,Ti for the total spin. In Table II, the product 
R\ i?2 is decomposed into the direct sum of irreducible representations. 

Table I. The number of each representation of 5*0(3, Z) which appears in the angular momentum 
L representation of SO(3,R). P = (— 1) L denotes the eigenvalue under parity transformation. 
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1 
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Table II. The decomposition of a product of two irreducible representations, Ri ® R 2 , into irre- 
ducible representations in 5*0(3, Z). Note that Ri ® R2 = R2 ® Ri by definition. 





A x 


A 2 


E 


Ti 






T 2 




Ai 


Ax 


A 2 


E 


Ti 






T 2 




A 2 


A 2 


Ax 


E 


T 2 






Ti 




E 
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E 


Ai ffi A 2 ®E 


Ti e t 2 






Ti ffiT 2 




Ti 


Ti 


T 2 


Ti ®T 2 


Ai ffi £ © Ti € 


DT 2 


A 2 ( 


B £ ffi Ti € 


E>T 2 


T 2 


T 2 


Tj 


Ti 0T 2 


A 2 ® £ ® Ti € 


BT 2 


Ai< 


B £ ffi Ti € 


DT 2 



We often use the wall source at t = to defined by 

J w&X \hUj 9 = BZf(t Q )Bjf{t Q ) (3-4) 

where a,/3 = 1,2 are upper component of the spinor indices while f,g are flavor 
indices. Here -B wall (tn) is obtained by replacing the local quark field q(x) of B(x) by 
the wall source, 

q ™n {to) = J2q(x,t ) (3-5) 

X 

with the Coulomb gauge fixing at t = to- Note that this gauge-dependence of the 
source operator disappears for the potential. All states created by the wall source 
have zero total momentum. Among them the state with zero relative momentum 
has the largest magnitude. A reason for employing the wall source here is that the 
ground state saturation for the potential at long distance is better achieved for the 
wall source than for other sources. 

Let us consider the case of the two nucleons. The source operator J waa {to) has 
zero orbital angular momentum at t = to, which corresponds to the A\ representation 
with positive parity. Therefore, the total angular momentum can be fixed by using 
the spin recoupling matrix M^ S ' Sz \ e.g., M^ s=0 ' Sz=0 "> = o- 2 /V2 and M^ s=1 ' Sz=m ^> = 
(o~20~ m ) I V% for m = 0, ±1 as 

J(t ; J P=+ , J z =m,I) = MfJj™ l \toU,f 9 . (3-6) 

Here P = ± is the parity and I = 1,0 is the total isospin of the system. Since the 
nucleon is a fermion, exchange of the nucleon operators in the source should give a 
minus sign. This fact fixes the total isospin given the total spin: (S, I) = (0, 1) or 
(1,0). (Note that S, I = are antisymmetric while S, I = 1 are symmetric under the 
exchange.) Since Af <g> A^S = 0) = Af and Af <g) Ti(S = 1) = T+, the state with 
either (J p , I) = (Af , 1) for the spin-singlet or (J p , I) = (T+,0) for the spin-triplet 
is created at t = to by the corresponding source operator. The NBS wave function 
extracted at t > to has the same quantum numbers (J P ,I) as they are conserved 
under QCD interactions. In addition the total spin S is conserved at t > to for the 
two nucleon system with equal up and down quark masses: Under the exchange of 
the two particles, the constraint (— l) s '+ 1 + / + 1 p = — 1 must be satisfied due to the 
fermionic nature of the nucleon. Also, the parity P and the isospin / are conserved 
in this system. Therefore S is conserved. However, L is not conserved in general. 
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While the state with (J p ,1) = (Af, 1) always has L = even at t > to, the one 
with (J P ,I) = (T+,0) has both L = A+ and L = E + ,T+ components*) at t > t , 
which corresponds to L = and L = 2 in 50(3, R), respectively. Note that J 
and L are used to represent the total and orbital quantum numbers respectively for 
50(3, Z) as well as for 50(3, R). 

The orbital angular momentum L of the NBS wave function for NN can be fixed 
to a particular value by the projection operator P^ as 

<p w {r; J p , I, L, 5) = P^p( s ^ w (r; J p , I) (3-7) 

where ip w (r; J p , I) is extracted from 

F(r, t - t ; J P ,I) ~ A{ J p , 1)^ \r; J p , I)e~ w ^\ (3-8) 
A{J P ,I) = (2B,W\J(t ;J p ,I)\0) 

for large t— to- The total spin projection operator is {P^^^a/B-^'p' = ^(c J 2)o/3(o'2)/3'o' 
for spin-sing let and p( s=1 ) = I — p( 5 =°) for spin-triplet, but this is redundant since 
the total spin 5, already fixed by the source, is conserved as mentioned before. The 
projection operator P( L ) of the orbital angular momentum for an arbitrary function 
ip(r) is defined in general by 

*% W lT) = ± £ X L {9)^ W (9- l -r) (3-9) 
g&SO(3,Z) 

for L = A\, A<i,E ,T\,T2, where x L denotes the character of the representation L in 
50(3, Z), * is its complex conjugate, g is one of 24 elements in 50(3, Z) and oIl is 
the dimension of L. 

3.2. Leading order NN potential: spin-singlet case 

We present the procedure to determine potentials at the leading order (LO): 

V LO (r) = V {r) + Kr(r)(<ri • <r 2 ) + V T (r)S 12 . (3-10) 

Since S±2 = and ct\ ■ a 2 = —3 for the spin-singlet case, the LO central potential 
for the spin-singlet case is extracted from the [J p ,1) = (Af, 1) state as 

Vr(r)^=W = yi=i (r) _ w i=i (r) _ [E k -H Q ]^(r;At,I = l,L = A 1 ,S = 0) 

(3-H) 

where 1 = + in isospin space. The potential Vc(r)( s,p) ^ 0,1 - ) in the above 
is often referred to as the central potential for the 1 So state, where the notation 
2S+1 Lj represents the orbital angular momentum L (see Table I), the total spin 5 
and the total angular momentum J of J = L + S. It is noted, however, that in 
the leading order of the velocity expansion, the potential does not depend on the 

*-* This can be seen from Table II for R2 = T\ (spin-triplet), which also tells us the existence of 
L — Tj + component in addition. The extra component is expected to be small since it appears as a 
consequence of the violation of 5*0(3, R) on the hyper-cubic lattice. 
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quantum number of the state J = L = A\. Moreover the A\ state may contain 
L = 4, 6, • • • components other than L = 0, though the L = component may 
dominate. Therefore it is more precise to refer to Vciv)^' 1 ^^' 1 ^ as the spin-singlet 
(isospin-triplet) central potential determined from the state with J = L = A\. 
A possible difference of spin-singlet central potentials between this determination 
and others such as the one determined from J = L = E gives an estimate for 
contributions from higher order terms in the velocity expansion. 

3.3. Leading order potential: spin-triplet case 

Both the tensor potential Vt and central potential Vq appear in the LO for the 
spin-triplet case. Let us consider the determination from (J P ,I) = (T+,0) state. 
The Schrodinger equation for this state becomes 

#o + V c (r f^ + V T (r)S 12 ] <p W (r; J P = T+,I) = E k <p w (r; J P = T+, /) (3-12) 

with (S, I) = (1,0), where the spin-triplet central potential is given by 

Vc (r)^ ee Vt°(r) + Vt°(r), Vt° = V$- 3VJ. (3-13) 

We separate the Schrodinger equation Eq. (3-12) into the A\ and non-^i components 
by using projection operators V = P^ 1 ) and Q = I — V as 

(Vc(r) {1 ' 0) - E k )V<p w {r) + V T (r)rS 12 <p w (r) = -H Vip w (r) 
(Vfc(r)(i.o) - E k )Q<p w (r) + V T {r)QS 12 <p w '(r) = -H Q<p w (r). (3-14) 

Note that V and Q commute with Hq, Vc(r) and Vr(r), whereas they do not com- 
mute with S\2- Non-^i component receives contributions from E, T\ and T 2 , among 
which only E and T 2 contribute to the D-wave. Since the contribution from T\ com- 
ponent turns out to be negligible in the numerical simulation, the non-Ai component 
is dominated by D-wave contributions. 

Using these projections, Vc and Vt can be extracted as 

V c (r) {1 ' 0) = E k-^r ) ([QSi2V W U(r)H [Pp W U(r) 

- [PS 12 ip w U(r)H [Q<p w U(r)) (3-15) 

VT{r) = ~W) (^ W U(r)H [Pp w U(r) - [V^ w U(r)H [Q^ w U(r))(3-W) 

A(r) ee [QS 12V w U(r)[V V w U(r) - [VS 12 <p w U(r)[Q<p w U(r). (3-17) 

In numerical simulations, (a, (3) = (2, 1) in J z = state is mainly employed. 

One may focus only on the A\ component of the wave function and define so- 
called the effective central potential for the spin-triplet (isospin-singlet), often used 
in nuclear physics: 
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The effect of Vr, which leads to a transition from the A\ component to the non- 
A\ component of the wave-function, is implicitly included in this effective central 
potential: For small Vr, the difference between Vc and Vq S is 0{Vj) as the second 
order perturbation tells us. 

3.4. Time- dependent HAL QCD method 

One of the practical difficulties to extract the NBS wave function and the poten- 
tial from the correlation function Eq.(3T) is to achieve the ground state saturation in 
numerical simulations at large but finite t— to with reasonably small statistical errors. 
While the stability of the potential against t — to has been confirmed within statis- 
tical errors in numerical simulations, 9 )' 18 ) the determination of W for the ground 
state suffers from systematic errors due to contaminations of possible excited states. 
There exist three different methods to determine W. The most well-known method 
is to determine W from the t — to dependence of the correlation function Eq.(3T) 
summed over r to pick up the zero momentum state. On the other hand, one may 
determine k 2 of W by fitting the r dependence of the NBS wave function with its 
expected asymptotic behavior at large r or by reading off the constant shift of the 
Laplacian part of the potential from zero at large r. Although the latter two methods 
usually give consistent results within statistical errors, the first method sometimes 
leads to a result different from those determined by the latter two at the value of 
t — to employed in numerical simulations. Although, in principle, the increase of 
t — to is needed in order to see an agreement among three methods, it is difficult in 
practice due to larger statistical errors at large t — to- 

The problem above is common in various applications of lattice QCD. Fortu- 
nately, the original HAL QCD method can be improved to overcome this difficulty 
as follows. Let us consider the normalized correlation function defined from Eq.(3-1) 
as 

= ^feS = Y,A n <p w »(r)exp(-tAW n ) + 0(e~ Aw ^), (3-19) 

n 

where AW n = W n — 2mB and AW^ = Wth — 2mB = tum- By neglecting the inelastic 
contributions above the meson production threshold, represented by 0(e~ AWtht ), for 
large enough t*\ non-relativistic approximation W n — 2m# ~ k^/niB leads us to 

R(r,t) ~ ^A n( p w "(r)exp (-*— ) = e~ t(Ho+u ^ A n ip w -{r) 

= e-^ Ho+u ^R(r,t = 0), (3-20) 

where the Schrodinger equation Eq. (2-5), the defining relation of the non-local 
potential U, is used to replace e~ tkn l niN by e _ *(^o+^)_ gy applying a time derivative 
on both side, we have the time-dependent Schodinger equation in imaginary time 

R(r, t)~J d 3 r'U(r, r')R(r', t). (3-21) 

*-* This limitation for t can be removed if the coupled channel potentials are introduced as in 
Sec. 7. 
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Now the velocity expansion of the non-local potential leads us to the formula of the 
leading order potential 

LO = (d/dt)R(r,t) _ H R(r,t) 

1 ' R{r,t) R(r,t) ' 1 ' 

Once the ground state saturation is achieved in R(r, t), Eq.(3-22) reduces to Eq.(2-17), 
for example for the spin-singlet case. Indeed, in this case, —d/dt is safely replaced 
by the non-relativistic energy of the ground state under the non-relativistic ap- 
proximation. 

The non-relativistic formula for V LO (r) above can be generalized to the case that 
masses of two particles are different by the replacement, R(r,t) = F(r,t)/e~( mA+mB ^ t . 
Note also that the potential extracted in this method automatically satisfies V LO (r — > 
) — > without constant shift. This property can be used to check whether this ex- 
traction works correctly or not. 

The non-relativistic approximation used to derive Eq.(3-19) can be removed by 
using the second order derivative in t; 



( 



\4m B dt 2 



^ - H ^j R(r, t) = j d 3 r'U(r, r')R(r', t), (3-23) 



which leads to 



V LO ( r ) - 1 (d/dtfR{r,t) _ (d/dt)R(r,t) _ H R(r,t) 

[ ' Am B R(r,t) R{r,t) R(r,t) ' 1 ' 

Here we have assumed that the inelastic contributions are negligibly small and that 
the two particles have the same mass. Since t is discrete on the lattice, the t deriva- 
tives has to be carefully performed. One has to employ the numerical derivative 
scheme which reduces statistical as well as systematic errors of V^°(r). 

One can generalize Eq.(3-19) to the correlation function with two relative coor- 
dinates x and y ; 

R ( x ,y,t) = —^i f d i x l d i y 1 {0\T{B(x l + x,t)B(x 1 ,t)B(y 1 + y,0)B( yi ,0)}|0), 

(3-25) 

which leads to 

4^W~di~ Ho ) R ^ y,t) = J d * z u ^ • z)jR(z ' y ' ^ (3 ' 26) 

Then, we obtain the non-local potential as 

d 3 z l-^g^ -m~ Ho ) R ^ z ^ • R~\ziVit)i (3-27) 
where R~ 1 (x,y,t) is an "truncated" inverse of the Hermitian operator R(x,y,t), 
R-\x,y,t) = ^j-rv n (x,t)vl(y,t) (3-28) 
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r [fm] r [fm] 



Fig. 2. (Left)The NN wave function for the spin-singlet and spin-triplet channels in the orbital ^4^~ 
representation at ~ 529 MeV and a ~ 0.137 fm in quenched QCD. The insert is a three- 
dimensional plot of the spin-singlet wave function ip w {x, y, z = 0). (Right) The NN (effective) 
central potential for the spin-singlet (spin-triplet) channel determined from the orbital Af wave 
function. Both figures are taken from Ref. 18). 

with X n (t) and v n (x,t) being the eigenvalues and corresponding eigenvectors of 
R(x,y,t), respectively. Note that zero eigenvalues are removed in the above sum- 
mation. Suppose we introduce a modified potential as 

U(x,y) = U(x,y) + ]T CnV n (x,t)v^(y,t). (3-29) 

A n =0 

Then it satisfies the same Schrodinger equation for all possible values of c n , the 
non-local potential is not unique as discussed before. 

§4. AT" AT" potential from lattice QCD 



4.1. Central potential in quenched QCD 

Let us first show results in the quenched QCD, where creations and annihilations 
of virtual quark-antiquark pairs are neglected: The standard plaquette gauge action 
is employed on a 32 4 lattice at the bare gauge coupling constant fi = 6/g 2 = 5.7. This 
corresponds to the lattice spacing a ~ 0.137 fm (1/a = 1.44(2) GeV), determined 
from the p meson mass in the chiral limit, and the physical size of the lattice L ~ 4.4 
fm. 9 ^ As for the quark action, the standard Wilson fermion action is used at three 
different values of the quark mass corresponding to the pion mass m n ~ 731, 529, 380 
MeV and the nucleon mass rn/v — 1560, 1330, 1200 MeV, respectively. 

Fig. 2(Left) shows the NBS wave functions for the spin-singlet and the spin- 
triplet channels in the orbital A\ representation at m w ~ 529 MeV. These wave 
functions are normalized to be 1 at the largest spatial point r ~ 2.2 fm. The 
central potential in the spin-singlet channel and the effective central potential in 
the spin-triplet channel extracted from the wave functions at ~ 529 MeV are 
shown in Fig. 2(Right). These potentials reproduce the qualitative features of the 
phenomenological iViV potentials, namely the repulsive core at short distance sur- 
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Fig. 3. The central potentials for the spin-singlet channel from the orbital Af representation at 
three different pion masses in quenched QCD. Taken from Ref. 18). 



rounded by the attractive well at medium and long distances. Prom this figure one 
observes that the interaction range of the potential is smaller than 1.5 fm, showing 
that the box size L ~ 4.4 fm is large enough for the potential. Labels 1 5o and 3 Si of 
the potentials in the figure represent the fact that potentials are determined from A± 
wave functions, which are dominated by the S-wave component. Note here that the 
lattice artifacts are expected to be large for potentials (as well as wave functions) at 
short distance such that r ~ 0(a). Therefore our results at short distance should be 
considered to be qualitative , not quantitative, and this caution should be applied 
to all of our results in this paper otherwise stated. The continuum extrapolation is 
necessary to predict short distance behaviors of potentials quantitatively. Indeed, 
BB potentials in the continuum limit are shown to diverge as r — > 0. 28 )' 29 ) 

In Fig. 3, NN central potentials in the spin-singlet channel are shown for three 
different pion masses. The repulsion at short distance and the attraction at medium 
distance are simultaneously enhanced as the pion mass decreases. 

4.2. Tensor potential in quenched QCD 

In Fig. 4(Left), we show the A\ and non-^4i components of the NBS wave func- 
tion obtained from the J p = (and J z = S z = 0) states at ~ 529 MeV. The 
non-j4i wave function is multivalued as a function of r due to its angular depen- 
dence. For example, the (a, j3) = (2, 1) component of the L = 2 part of the non-Ai 
wave function is proportional to the spherical harmonics Y2o(6,4>) oc 3 cos 2 9 — 1. 
Fig. 4(Right) shows the non-yli component divided by Y2o(6,4>). The non-^i wave 
function seems to be dominated by the D (L = 2) state, since its multivaluedness 
is mostly absorbed to Y2q(9,4>). Fig. 5 (Left) shows the central potential Vc(r)( 1,0 ^ 
and tensor potential Vr(r), together with the effective central potential V^ ff (r)( 1,0 ), 
at the leading order of the velocity expansion as given in Eqs. (3T5), (3-16) and 
(3T8), respectively. 

Note that V(^(r) contains the effect of Vt(t) implicitly as higher order effects 
through processes such as 3 S\ — > 3 D\ — > 3 S\. At the physical pion mass, V^ ff (r) is 
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Fig. 4. (Left) (a, /3) = (2, 1) components of the orbital A{ and non-A^ wave functions from J p — 
7\ + (and J z — S z = 0) states at m n — 529 MeV. (Right) The same wave functions but the 
spherical harmonics components are removed from the non- A+ part. Taken from Ref. 18). 



100 




Fig. 5. (Left) The central potential Vc(?") °' and the tensor potential Vt(t) obtained from the 
J p = T+ NBS wave function, together with the effective central potential V§ a (rp 1 ' ', at 
— 529 MeV. (Right) Pion mass dependence of the tensor potential. The lines are the four- 
parameter fit using one-pion-exchange + one-rho-exchange with Gaussian form factor. Taken 
from Ref. 18). 



expected to gain sufficient attraction from the tensor potential, which leads to the 
appearance of a bound deuteron in the spin-triplet (and flavor-singlet) channel while 
an absence of the bound dineutron in the spin-singlet (and flavor-triplet) channel. 
The difference between Vc{r)^^ and V^ ff (r) in Fig. 5 (Left) is still small in this 
quenched simulation due to relatively large pion mass. 

The tensor potential in Fig. 5 (Left) is negative for the whole range of r within 
statistical errors and has a minimum around 0.4 fm. If the tensor potential receives a 
significant contribution from one-pion exchange as expected from the meson theory, 
Vt(t) is rather sensitive to the change of the pion mass. As shown in Fig. 5 (Right), 
it is indeed the case: Attraction of Vr(r) is substantially enhanced as the pion mass 
decreases. 

The central and tensor potentials in lattice QCD are given at discrete data 
points. For practical applications to nuclear physics, however, it is more convenient 
to parameterize the lattice results by known functions. Such a fit for Vt(?) is given 



Lattice QCD approach to Nuclear Physics 



17 



600 



° 200 



50 


45[MeV] — 
0[MeV] — . 


1 











^■^^ 




-50 

V 






1 2 





1 

r[fm] 



> 400 



Vc(r): D 2 (S=0 L=2) 
Vc(r): 1 S o (S=0 L=0) 

100 n . . .- 



- 




0.5 1 1.5 2 



0.5 



1 

[fm] 



1.5 



Fig. 6. (Left) The spin-singlet central potential VcM^ obtained from the orbital channel 
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Fig. 7. (Left) The spin-triplet central potential Voir) obtained from the orbital A\ — T 2 + cou- 
pled channel in quenched QCD at ~ 529 MeV. (Right) The tensor potential Vr(r) from the 
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For these two figures, symbols are same as in Fig. 6(Left). 



by the form of one-pion-exchange + one-rho-exchange with Gaussian form factors as 



V T (r) = bl (l-e-^y (1 + — + ^) 

\ m p r [m p r) A J 



m p r (m p r) 2 J r 

1 hi 1 - < -^ 2 ) 2 ( 1 + — + r^—] — , (4-i) 

m^r [m^ry J r 



where &i, 2,3,4 ar e the fitting parameters while m n (m„) is taken to be the pion mass 
(the rho meson mass) calculated at each pion mass. The fit line for each pion 
mass is drawn in Fig. 5 (Right). It may be worth mentioning that the pion-nucleon 
coupling constant extracted from the parameter 63 in the case of the lightest pion 
mass (mj,- = 380 MeV) gives g^ N /(A-K) = 12.1(2.7), which is encouragingly close to 
the empirical value. 
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4.3. Validity of velocity expansion in quenched QCD 

The potentials so far are derived in the leading order of the velocity expansion. 
It is therefore important to investigate the convergence of the velocity expansion. If 
the non-locality of the NN potentials were absent, the leading order approximation 
for the potentials would give exact results at all energies below the inelastic thresh- 
old. The non-locality of the potentials therefore becomes manifest in the energy 
dependence of the potentials. 

To study the energy dependence, the leading order local potentials at E ~ 45 
MeV, realized by anti-periodic boundary conditions in the spatial directions, are 
calculated in quenched QCD at m n ~ 529 MeV and L ~ 4.4 fm. 30 ^ 33 ^ In this 
calculation, four types of momentum-wall sources, defined by 

?(to;/) = 2>(x,t )/(x) (4-2) 

x 

are employed, where /(x) = cos((±x ± y + z)ir/L). Note that /(x) = 1 corresponds 
to the wall source used in the periodic boundary condition. These momentum-wall 
sources induce L = as well as L = A\ states. 

In Fig. 6(Left), the spin-singlet potential Vcir)^' 1 ^^ ' 1 ^ obtained from the 
L = A^ state at E ~ 45 MeV (red circles) is compared with that at E ~ MeV 
(blue circles), while a comparison is made in Fig. 7 for the spin-triplet potentials, 
V r c .(r)( 5 ' / ) = ( 1 '°^(left) and Vr(r) (right). Good agreements between results at two 
energies indicate that higher order contributions are rather small in this energy in- 
terval. In other words, these local potentials obtained at E ~ MeV can be safely 
used to describe the NN scattering phase shift in the range between E = MeV 
and E = 45 MeV at this pion mass in quenched QCD. 

Non-locality of the potential may become manifest also in its angular momentum 
dependence, since the orbital angular momentum L = r xp contains a derivative. In 
Fig. 6 (Right), the spin-singlet potential Vcir)^' 1 ^^ ' 1 ^ obtained from the L = 
state, whose main component has L = 2, is compared to the one from the L = 
Af state, whose main component has L = 0. In this comparison, local potentials 
are determined at the same energy, E ~ 45 MeV, but different orbital angular 
momentum. Although the statistical errors are rather large in the case of L = T^~, a 
good agreement between the two is again observed, suggesting that the L dependence 
of the potential is small at least for the spin-singlet case. 

By these comparisons, it is observed that both energy and orbital angular mo- 
mentum dependencies for local potentials are very weak within statistical errors. We 
therefore conclude that contributions from higher order terms in the velocity expan- 
sion are small and that the LO local potentials in the expansion obtained at E ~ 
MeV and L = are good approximations for the non-local potentials at least up to 
the energy E ~ 45 MeV and orbital angular momentum L = 2. 

4.4. Central potential in full QCD 

Needless to say, it is important to carry out calculations of iViV potentials in 
full QCD on larger volumes at lighter pion masses. The PACS-CS collaboration is 
performing 2 + 1 flavor QCD simulations, which cover the physical pion mass. 34 )' 35 ) 
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Fig. 8. (left) The multi-Gaussian fit of the central potential Vc(r) with iVcauss = 5. (right) The 
scattering phase in 1 5'o channel in the laboratory frame obtained from the lattice NN potential, 
together with experimental data 38). 



Gauge configurations are generated with the Iwasaki gauge action and 0(a)-improved 
Wilson quark action on a 32 3 x 64 lattice. The lattice spacing a is determined from 
m n , ttik and tjiq as a ~ 0.089 fm, leading to L ~ 2.9 fm. Three ensembles of 
gauge configurations are used to calculate NN potentials at (m^rn/v) —(701 MeV, 
1583 MeV), (570 MeV, 1412 MeV) and (411 MeV,1215 MeV ) 36 > . To overcome 
a difficulty to achieve ground state saturations in full QCD simulations, the time- 
dependent HAL QCD method in Sec. 3.4 is employed. 37 ^ 

Fig. 8(Left) shows the spin-singlet NN central potential Vc(r) obtained at 
E ~ from the PACS-CS configurations with ~ 701 MeV and — 1583 
MeV. This central potential Vc(r) is fitted with multi-Gaussian function that g(r) = 
^2n=i uss ex P( — v n r2 ) with fit parameters V n and u n {> 0). A solid line in the figure 
represents a fit result with A^uss = 5. 

We solve the Schrodinger equation in 1 Sq channel with this fitted potential Vc(r), 
in order to calculate the scattering phase shift. Fig. 8(Right) shows the scattering 
phase S(k) in the laboratory frame, together with the experimental data 38 ) for a 
comparison. A qualitative feature of the experimental data is well reproduced by 
the lattice potential, though the strength is weaker, most likely due to the heavier 
pion mass, m n ~ 701 MeV. The scattering length obtained from the derivative of 
the phase shift at k = becomes a( 1 5o) = lim/^o tan 5(k)/k = 1.6(1.1) fm, which is 
compared to the experimental value a exp ( 1 S'o) — 20 fm. 

4.5. Nuclear force in odd parity sector and the spin-orbit force in full QCD 

In this subsection, we consider the potentials in odd parity sectors. Together 
with the nuclear forces in even parity sectors, the information on odd parity sectors 
is necessary in studying many-nucleon systems with Schrodinger equations. In par- 
ticular, we are interested in the spin-orbit (LS) force, which gives rise to a part of 
the spin-orbit coupling in the average single-particle potential of nuclei. It is also 
expected to induce superfluidity in neutron stars by providing an attraction between 
two neutrons in 3 P2 channel. 13 ^ 
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The LS force appear at NLO of the derivative expansion as 

#o + V c (r)^ I} + V T (r)S 12 + V LS (r)L ■ S 

To obtain three unknown potentials, Vc, Vr and Vls, we need three independent 
NBS wave functions. We therefore generalize the two-nucleon source for odd parity 
sectors, by imposing a momentum on the composite nucleon fields as 

J a p(t ; f^) ^ N a (t ; f^)Np(t ; f^*) for i = ±l,±2,±3, (44) 

where TV denotes a composite nucleon source field carrying a momentum, 

N a (t ; /«) = Y, tabc (uI(x 1 )C 75 4(x 2 )) g c , a (x 3 )/»(x 3 ), (4-5) 

Xl,X 2 ,X 3 

with /( ±3 ')(x) = exp[±27TKCj/L]. The star "*" in the r.h.s. of Eq. (4-4) represents 
the complex conjugation, which is used to invert the direction of the plane wave. A 
cubic group analysis shows that the two-nucleon source Eq. (4-4) contains the orbital 
contribution Af © E + © XT/, whose main components are S-wave, D-wave and P- 
wave, respectively. Thus the two-nucleon source Eq. (4-4) covers all the two-nucleon 
channels with J < 2. 

For the spin-triplet odd-parity sector, Eq. (4-4) generates the lowest-lying NBS 
wave functions for (J P ,I) = (A± , 1), (T-j - , 1), (E~ , 1) and (T 2 ~,l), which roughly 
correspond to J p = 0~,1~,2~ and 2~, respectively. Among these, we consider 
Schrodinger equations for three NBS wave functions in J p = A± ,T± ,T 2 ~ as 

[H + V c (r) + V T (r)S l2 + V LS (r)L ■ S] ip w (r; J p ) = E (J p )ip w (r; J p ), (4-6) 

where Eq(J p ) = k 2 /mjv from the lowest-lying energy W = 2^J m? N + k 2 for the J p 
sector. In order to obtain Vc(r), Vr(r) and Vls{t) in odd parity sectors, Eqs. (4-6) 
for J p = A^,T{ ,T 2 , which correspond to 3 P , 3 Pi and 3 P 2 + 3 F 2 , are solved. 

Numerical calculations are performed by using 2 flavor QCD gauge configura- 
tions on 16 3 x 32 lattice generated by CP-PACS Collaboration, 39 )- 40 ) with Iwasaki 
gauge action at /3 = 1.95 and 0(a) improved Wilson (clover) quark action at 
k = 0.1375. This setup leads to the lattice spacing a^ 1 = 1.27 GeV (a ~ 1.555 
fm), the pion mass ~ 1136 MeV, the nucleon mass ~ 2165 MeV. The spatial 
extension amounts to L = 16a ~ 2.5 fm. 

Fig. 9 shows preliminary results of the central potential Vc(r), tensor potential 
Vr(r) and the spin-orbit force Vls{t) in the spin-triplet odd parity sector. They have 
the following qualitative features. (1) Vc(r) has repulsive core at short distance. (2) 
Vr(r) is positive and very small. (3) Vls(t) is large and negative at short distance. 
These features qualitatively agree with those of phenomenological potentials. 8 ) 

§5. Hyperon Interactions 

Study of hyperon-nucleon (YN) and hyperon-hyperon (YY) interactions is one 
of the challenges in contemporary nuclear physics. These potentials give a key to 



<p w (r;J-,I) = E k <p w {r;J-,I) (4-3) 
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Fig. 9. Potentials in odd parity sector obtained from 3 Po, 3 Pi and 3 P2 + 3 F2 NBS wave functions 
calculated at — 1136 MeV. Left, middle, and right figures show central, tensor and spin-orbit 
force in parity odd sector, respectively. 



understand nuclear many-body systems with strangeness. Also they are essential to 
explore the structure of the neutron star core, where strangeness degree of freedom 
is expected to appear. At present, experimental data on YN and YY scatterings 
are not sufficient to make precise constraints on the hyperon potentials, while spec- 
troscopic studies of A hypernuclei, performed by employing various reactions such 
as (tt + ,K + ), (K~,tt~) and (e, e'X + ), 41 ^ 42 ^ give some information on the AN in- 
teractions. Under these circumstances, studies on the basis of lattice QCD is quite 
important as an alternative method to access YN and YY interactions. In this sec- 
tion we mainly consider potentials in the strangeness S = — 1 sector, obtained from 
2 + 1 flavor lattice QCD simulations with PACS-CS gauge configurations. A study 
on potentials between octet barons in the flavor SU(3) limit and coupled channel 
analysis on potentials in the strangeness S = —2 sector beyond the SU(3) limit will 
be discussed in the next two sections. 

5.1. AN and UN potentials in full QCD 

The AN and the UN (I = 3/2) are the lowest states in the strangeness S = — 1 
systems with 1 = 1/2 and / = 3/2, respectively. Therefore potentials for these states 
can be calculated as in the case of NN potential. In Ref. 43), the AN potential 
and the UN potential with / = 3/2 are calculated by using 2+1 flavor full QCD 
gauge configurations with the original time-independent HAL QCD method. In the 
following, we show improved results on a 32 3 x 64 lattice at a = 0.091(1) fm with 
the time-dependent HAL QCD method discussed in Sec. 3.4. 

The AN (left panel) and the UN(I = 3/2) (right panel) potentials in the 1 S 
channel are shown in Fig. 10. In the 2+1 flavor QCD, while the UN (I = 3/2) poten- 
tial still belongs directly to the 27(7 = 3/2) representation thanks to the isospin sym- 
metry, an energy eigenstate of a AN system in the 1 5o channel is a mixture of 27(1 = 
1/2) and 8 S in the flavor representation, so that these two potentials are not necessar- 
ily equal. In the present 2 + 1 flavor QCD calculation shown in Fig. 10, these poten- 
tials look similar due to small flavor-SU(3) breaking: For example, our hadron masses 
are (m n ,m K ,m N ,m A ,m E ) = (0.7006(4), 0.7879(4), 1.574(3), 1.635(3), 1.650(3)) GeV. 
The left panel of Fig. 11 shows the central potential (circle) and the tensor potential 
(triangle) of the AN system in the 3 Si — 3 D\ channel, whose eigenstate is a mixture 
of 10 and 8 a . The attractive well at distance r ~ 0.6 fm is deeper than that of the 
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Fig. 10. Left: The central potential in the So channel of the AN system in 2 + 1 flavor QCD as a 
function of r. Right: The central potential in the 1 So channel of the UN (I = 3/2) system as a 
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AN central potential in the 1 5'o channel, while the tensor potential itself (triangle) 
is weaker than the tensor potential in the NN system. 44 ) 

The right panel of Fig. 11 shows the central potential (circle) and the tensor 
potential (triangle) of the EN (I = 3/2) system in the 3 5i — 3 D\ channel. Due 
to the isospin symmetry, this channel belongs solely to the flavor 10 representation 
without mixture of 10 or 8 a As seen from the figure, there is no clear attractive well in 
the central potential (circle). This repulsive nature of the EN (I = 3/2, 3 Si — 3 D\) 
central potential is consistent with the prediction from the naive quark model. 45 ) 
The tensor force is a little stronger that that of the AN system but is still weaker in 
magnitude than that of the NN system. 

5.2. SN potential in quenched QCD 

Experimentally, not much information is available on the NS interaction ex- 
cept for a few studies: a recent report gives the upper limit of elastic and inelastic 
cross sections 46 ^ while earlier publications suggest weak attractions of S— nuclear 
interactions. 47 )~ 49 ) The S— nucleus interactions will be soon studied as one of the 
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Fig. 12. (Left) The spin-singlet central potential for pS° obtained from the orbital channel at 
~ 368 MeV (circle) and m„ ~ 511 MeV (box). The central part of the OPEP (F/(F + D) = 
0.36) in Eq. (5-2) is also given by solid line. (Right) The spin-triplet effective central potential 
from the orbital Af channel at m n — 368 MeV (triangle) and m n — 511 MeV (diamond), 
together with the OPEP (solid line). Taken from Ref. 51). 



day-one experiments at J-PARC 50 - 1 via (K~,K + ) reaction with a nuclear target. 
Ref. 51) gives the first result of the potential in I = 1 NS system, which does not 
show strong decay into other channels. Lattice parameters are the same as for the 
quenched NN potential in Sec.4.1, but the method to determine the lattice spacing 
in Ref. 51) is a little different from the one in Sec.4.1. The potential is calculated 
at (m^,m N ,m s ) = (511(1) MeV, 1300(4) MeV, 1419(4) MeV) and (368(1) MeV, 
1167(7) MeV, 1383(6) MeV) with the interpolation operators 

p a (x) = e abc (ul(x)C-f 5 d b (x))u Cta (x), = e ahc {u^{x)C^s h {x))s Cta {x). (5-1) 

Since both p and S° have (I,I Z ) = (1/2, 1/2), the pS° system has 1 = 1 with the 
strangeness S = —2. 

The left (right) of Fig. 12 gives the (effective) central potential of the pS° system 
obtained from the L = Af representation for the spin-singlet (triplet) at m n = 511 
MeV and 368 MeV. Potentials in the 1 = 1 NS system for both channels show a 
repulsive core at r < 0.5 fm surrounded by an attractive well, similar to the NN 
systems. In contrast to the NN case, however, the repulsive core of the pS° potential 
in the spin-singlet channel is substantially stronger than in the triplet channel. The 
attraction in the medium to long distance region ( 0.6 fm < r < 1.2 fm ) is similar 
in both channels. The height of the repulsive core increases as the light quark mass 
decreases, while a significant difference is not seen for the attraction in the medium 
to long distance within statistical errors. Potentials in Fig. 12 are weakly attractive 
on the whole in both spin channels at both pion masses, in spite of the repulsive 
core at short distance, and the attraction in the triplet is a little stronger than that 
in the singlet. 

The solid lines in Fig. 12 are the one-pion exchange potential (OPEP), given by 

VS = -(1 - 2a) ^ (t n .t 5 )(«n-<ts) fm^Ve^ 
c y 1 4vr 3 \2m N J r v 1 

with (m,r, mjv) = (368MeV, 1167MeV), where the pseudo- vector ttSS coupling f^ss 
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Fig. 13. The BB potentials in 27 (Left) and 10 (Right) representations extracted from the lattice 
QCD simulation at M ps = 469 MeV. Taken from Ref. 56). 



is related to the irNN coupling as f n ~~ = ~/irJVJv(l — 2a) with the parameter a = 
F/(F + D), and g n NN = fnNN^^-- The empirical vales, a ~ 0.36 and g n NN/(^) — 
14.0, are used for the plot. Unlike the NN potential, the OPEP in the present 
case has opposite sign between the spin-singlet channel and spin-triplet channel. 
In addition, the absolute magnitude is smaller due to the factor 1 — 2a. No clear 
signature of the OPEP at long distance (r > 1.2 fm) is yet observed in Fig. 12 within 
statistical errors. 



§6. Baryon interaction in the flavor SU(3) limit 



6.1. Potentials in the flavor SU(3) limit 

In order to reveal the nature of the hyperon interactions in various channels, 
it is more convenient to consider an idealized flavor SU(3) symmetric world, where 
u, d and s quarks are all degenerate with a common finite mass. In this limit, one 
can capture essential features of the interaction, in particular, the short range force 
without contamination from the quark mass difference. 

In the flavor SU(3) limit, the ground state baryon belongs to the flavor-octet 
with spin 1/2, and two-baryon states with a given angular momentum can be labeled 
by the irreducible representation of SU(3) as 

8<g>8 = 27e8eieIoeioe8, (6-1) 

V v ' S v ' 

symmetric anti— symmetric 

where "symmetric" and "anti-symmetric" stand for the symmetry under the ex- 
change of the flavor for two baryons. For the system with orbital S-wave, the Pauli 
principle for baryons imposes 27, 8 and 1 to be spin-sing let (^o), while 10, 10 and 
8 to be spin-triplet ( 3 Si — 3 D\). Calculations in the SU(3) limit allow us to extract 
potentials for these six flavor irreducible multiplets as follows. 

A two-baryon operator BB^ which belongs to one definite flavor representation 
X, can be given in terms of the baryon base operator with the corresponding Clebsch- 
Gordan (CG) coefficients C% as BB^ = Y.ij C vj B i B r B Y 

using this operator at 

source and/or sink, the NBS wave function for two-baryon system in the flavor 



Lattice QCD approach to Nuclear Physics 



25 



1500 
1000 
500 




S.v) 



6000 
5000 
4000 
3000 
2000 
1000 




Vc 



>c ua . s =0. 13840 



0.0 0.5 1.0 1.5 2.0 2.5 



1500 
1000 
500 



10) 




K, a 5=0.13840 



0.0 0.5 1.0 1.5 2.0 2.5 



0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 



yd) 























i 










-50 










-100 










-150 




ic Ui a lS =0. 13840 






0.0 0.5 1.0 


1.5 2.0 2.5 



8a) 



Vi > 




-100 - Kji,a.s=0.1384(l . 

0.0 0.5 1.0 1.5 2.0 2.5 



0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 

r [fm] r [fm] 
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representation can be obtained. Potentials in the flavor base, K"(r), are extracted 
form such wave functions in the same manner for nuclear forces explained in Sec. 3. 

In Ref.52), the (effective) central potentials are calculated in the original time- 
independent HAL QCD method by using the 3 flavor full QCD gauge configuration 53 ) 
on a 16 3 x 32 lattice at a ~ 0.12 fm, and at two values of quark hopping parameter cor- 
responding to (M ps ,M B ) = (1014(l)MeV,2026(3)MeV) and (835(1) MeV, 1752(3) 
MeV), where M ps and Mb denote the octet pseudo-scalar (PS) meson mass and the 
octet baryon mass, respectively. In Refs. 54)— 56), on the other hand, the central and 
tensor potentials are calculated in the time-dependent HAL QCD method discussed 
in Sec. 3. 4 by using the 3 flavor full QCD gauge configuration generated by HAL 
QCD Collaboration on a 32 3 x 32 lattice at a ~ 0.12 fm, and at five values of quark 
hopping parameter which correspond to (M ps , M B ) = (1170.9(7)MeV, 2274(2)MeV), 
(1015(1) MeV, 2030(2) MeV ), (837(1) MeV, 1748(1) MeV ), (673(1) MeV, 1485(2) 
MeV ), and (468.6(7) MeV, 1161(2) MeV ). Figs. 13 and 14 show the flavor basis 
potentials for M ps = 469 MeV. 56 ) The left panels show central potentials for the 
spin-singlet channel from the J p = Af state, while the right panels give central (Vc) 
and tensor (Vr) potentials for the spin-triplet channel from the J p = T x + state. 

As listed in Table III, some of octet-baryon pairs solely belong to an irreducible 
representation of flavor SU(3). For example, symmetric AW belongs to 27 represen- 
tation. Therefore, V^ 27 \r) can be considered as flavor SU(3) symmetric limit of the 
NN spin-sing let ^So) potential. Similarly V^, V^ and V( 8 -> can be considered 
as flavor SU(3) symmetric limit of some BB potentials of the particle basis, while 
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flavor multiplet 


baryon pair (isospin) 


27 


{JVJV}(I=1), {iVr} (1=3/2), {EE}(1=2), 




{JC~}(I=3/2), {SS}(I=1) 


8 S 


none 


1 


none 


10* 


[NN] (1=0), [JC~] (1=3/2) 


10 


[NS] (1=3/2), [SS](I=0) 


8a 


[ATS] (1=0) 



Table III. Baryon pairs which belongs to an irreducible flavor SU(3) representation, where {BB 1 } 
and [BB'] denotes BB' + B'B and BB' - B'B, respectively. 

a nd are always mixtures of different BB potentials in the particle basis. 

Fig. 13 shows V^? 7 \r) and VJ^P (r), which correspond to spin-singlet and spin- 
triplet NN potentials, respectively. Both central potentials have a repulsive core at 
short distance with an attractive pocket around 0.8 fm. These qualitative features 
are consistent with the results found for the NN potential in previous section. The 
upper-right panel of Fig. 14 shows that V^. 10 ^ (r) has a stronger repulsive core and a 

weaker attractive pocket compared to Vq (r) and Vq (r). Furthermore Vq (r) 
in the upper-left panel of Fig. 14 has a very strong repulsive core among all 6 channels, 

(8 ) 

while V(j (r) in the lower-right panel has a very weak repulsive core. In contrast 
to all other cases, V^\r) has attraction at short distances instead of repulsion, as 
shown in the lower-left panel. 

Above features are consistent with what has been observed in a SU(6) quark 
model. 45 ) In particular, the potential in the 8 S channel in this quark model becomes 
strongly repulsive at short distance since the six quarks cannot occupy the same 
orbital state due to the Pauli exclusion for quarks. On the other hand, the potential 
in the 1 channel does not suffer from the quark Pauli exclusion at all, and can 
become attractive due to the short-range gluon exchange. Such agreements between 
the lattice data and the quark model suggest that the quark Pauli exclusion plays 
an essential role for the repulsive core in BB systems. 

The potential for the flavor singlet is entirely attractive even at very short dis- 
tance. This might produce a bound state, the //-dibaryon, in this channel. We will 
discuss this possibility in the next subsection. 

In the flavor SU(3) limit, the BB potentials in the particle basis can be obtained 
from those in flavor basis by a unitary rotation as 

Viiir) = Y J U ixV^ X \r)U* Xj (6-2) 
x 

where Uix is a unitary matrix which rotates the flavor basis \X) to the particle basis 
\i) as \i) = Uix\X), and given in terms of the CG coefficients. The explicit forms of 
the unitary matrix U are found in Ref. 52). 

In Fig. 15, we show BB potentials for S=-2, 1=0, sector at M PS = 469 
MeV, as a characteristic example. The flavor base potentials are fitted by the ana- 
lytic function composed of an attractive Gaussian core plus a long range (Yukawa) 2 
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(Center) and spin-triplet Vt (Right), extracted from the lattice QCD simulation at M ps = 469 
MeV. 



attraction, 

with five parameters ^1,2,3,4,5- The left panel of Fig. 15 shows the diagonal potentials. 
One observes that all three diagonal potentials have a repulsive core. The repulsion 
is most strong in the 1717(1=0) channel, reflecting its largest CG coefficient of the 8 S 
state among three channels, while the attraction in the 1 state is reflected most in 
the ^^(1=0) potential due to its largest CG coefficient. The right panel of Fig. 15 
shows the off-diagonal potentials, which are comparable in magnitude to the diagonal 
ones, except for the AA-NS transition potential. Since the off-diagonal parts are 
not negligible in the particle basis, a fully coupled channel analysis is necessary to 
study observables in this system. This is important when we study the real world 
with the flavor SU(3) breaking, where only the particle base is meaningful. 

As another example, we show BB potentials for S=— 1, 1=1/2 sector at Mpg = 
469 MeV in Fig. 16. The left panel of Fig. 16 shows the potential in the 1 Sq channel, 
while the center (right) panel of Fig. 16 shows the central (tensor) potential for 
J p = T+ spin-triplet channel. We observe that the off-diagonal NA-NE potentials 
are significantly large, especially in the spin-triplet channel, so that the full NA-NE 
coupled channel analysis is also necessary to study observables. In addition, the 
repulsive cores in the spin-single channel are much stronger than that in the spin- 
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triplet channel, due to the strong repulsion of the 8 S state, which couples only to 
the spin-singlet channel. 

Although all quark masses of 3 flavors are degenerate and rather heavy in the 
present simulations, these particle base potentials, shown in Fig. 15 and Fig. 16, 
may provide useful information for the behavior of hyperons in hyper-nuclei and in 
neutron stars. 42 )' 57 ) BB potentials in all other sectors can be found in Ref. 52). 

6.2. Bound H dibaryon in the flavor SU(3) limit 

In this subsection, we investigate the potential of the flavor singlet channel in 
order to see whether the bound H dibaryon exists or not in the flavor SU(3) limit 
case. 

Potentials for the flavor irreducible channels in the SU(3) limit have been calcu- 
lated in Ref. 54)-56) on 16 3 x 32, 24 3 x 32 and 32 3 x 32 lattices at a = 0.121(2) fm 
and five values of the quark mass, as mentioned before. 

Shown in Fig. 17(Left) and Fig. 17(Right) are the volume and the quark mass 
dependencies of the central potential in the flavor-singlet channel Vq(t) at (t — 
to) /a = 10, where the potentials do not have appreciable change with respect to 
the choice of t. The flavor-singlet potential is shown to have an "attractive core" 
and to be well localized in space. Because of the latter property, no significant 
volume dependence of the potential is observed within the statistical errors, as seen 
in Fig. 17(Left). As the quark mass decreases in Fig. 17(Right), the long range part 
of the attraction tends to increase. 

The resultant potential is fitted by the form in Eq. (6-3). With the five parame- 
ters, 61,2,3,4,5, the lattice results can be fitted reasonably well with x 2 /dof ~ 1. The 
fitted result for L = 3.87 fm is shown by the dashed line in Fig. 17(Left). 

Solving the Schrodinger equation with the fitted potential in infinite volume, 
the energies and the wave functions are obtained at the present quark masses in the 
flavor SU(3) limit. It turns out that, at each quark mass, there is only one bound 
state with binding energy of 20-50 MeV. Fig. 18(Left) shows the energy and the 
root-mean-squared (rms) distance of the bound state at each quark mass obtained 
from the potential at L = 3.87 fm and (t — to)/a = 10, where errors are estimated by 
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the jackknife method. Despite the fact that the potential becomes more attractive 
as quark mass decrease, the resultant binding energies of the i/-dibaryon decrease 
in the present range of the quark masses, since the increase of the attraction toward 
the lighter quark mass is compensated by the increase of the kinetic energy for the 
lighter baryon mass. It is noted that there appears no bound state for the potential 
of the 27-plet channel or the 10-plet channel ("deuteron" ) in the present range of 
the quark masses. 

By including systematic errors caused by the choice of sink-time t in R(r, t — to), 
the final results of the binding energy Bh and the rms distance \J (r 2 ) are summarized 
below, where the 1st and 2nd parentheses correspond to statistical and systematic 
errors, respectively. 



M ps = 


1171 MeV : 


B H 


= 49.1(3.4)(5.5) MeV y/(r 2 ) 


= 0.685(13)(25) fm 


M ps = 


1015 MeV : 


B h 


= 37.2(3.7)(2.4) MeV V(r 2 ) 


= 0.809(23)(10) fm 


M ps = 


837 MeV : 


Bh 


= 37.8(3.1)(4.2) MeV vV> 


= 0.865(20)(25) fm 


M ps = 


672 MeV : 


Bh 


= 33.6(4.8)(3.5) MeV ^/(r 2 ) 


= 1.029(41)(23) fm 


M ps = 


469 MeV : 


B h 


= 26.0(4.4)(4.8) MeV ^(r 2 ) 


= 1.247(70)(59) fm. 



Recently, the existence of iJ-dibaryon is also investigated by a direct calculation 
of its binding energy in 2+1 full QCD simulations, 58 )' 59 ) where B H = 13.2(1.8)(4.0) 
MeV is reported in the L — > oo extrapolation at ~ 389 MeV, wik — 544 MeV. 
Fig. 18(Right) gives a summary of the .ff-dibaryon binding energy from full QCD 
simulations recently reported. 

Since the binding energy is comparable to the splitting between physical hyperon 
masses and not so sensitive to quark mass, there may be a possibility of weakly 
bound or resonant ii-dibaryon even in the real world with lighter quark masses and 
the flavor SU(3) breaking. Our phenomenological trial analysis using 3-flavor lattice 
QCD results, suggests a resonant If-dibaryon above AA but bellow NS thresholds. 56 ) 
To make a definite conclusion on this point, however, the AA — NS — EE coupled 
channel analysis is necessary for H in the (2+l)-flavor lattice QCD simulations, as 
will be discussed in Sec. 7. 

§7. Hadronic interactions above inelastic threshold 

In this section, we discuss a method to investigate hadron interactions above 
inelastic threshold by generalizing the HAL QCD method. Then, we apply the 
method to coupled channel potentials in the S = — 2 and 1 = sector. 

7.1. Coupled- channel approach to inelastic scattering 

For simplicity, let us discuss a case oi A+B ^ C+D scattering where A,B,C,D 
represent some 1-particle states. This is a simplified version of the octet baryon 
scattering in the strangeness S = — 2 and isospin 1 = channel, where AA, NS and 
EE appear as asymptotic states of the strong interaction if the total energy is larger 
than 2ms- We here assume tua + mg < mc + mo < W, where W = + 
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Fig. 18. (Left) The energy Eo and the root-mean-square distance 



2 ) of the bound state in 
the flavor singlet channel at each quark mass. Bars represent statistical errors only. (Right) 
Summary of the Jf-dibaryon binding energy in recent full QCD simulations. HAL stands for 
the present results and NPL stands for the result in Ref. 59 ' 



is the total energy of the system, and = yra| + k 2 . In this case, the QCD 
eigenstate with the quantum numbers of the AB state and center of mass energy W 
is expressed as 

\W) = c AB \AB, W) + c CD \CD, W) + --- (7-1) 
\AB,W) = \A,k) in ®\B,-k)- m , \CD,W) = \C,q) in ®\D,-q} in , (7-2) 

where W = E£ + Eg = E° + Then we define the following NBS wave functions, 

tpABir, k)e~ wt = (0\T{ip A (x + r, t)cp B (x, t)}\W), (7-3) 
fCD(r, q)e~ wt = (0|T{ W (x + r, t)<p D (x, t)}\W). (7-4) 

Using the partial wave decomposition such that*) 

<px(r,k) = 4tt^2 i e ip i x (.r,k)Y em (f2 r )Y im {f2 k ) (7-5) 

l,m 

for X = AB or CD, it can be shown 61 ) that these wave functions satisfy 

{V 2 + k 2 ) VAB {r,k) =0, (V 2 + q 2 )<p CD (r,q) =0 (7-6) 

for r — > oo. 

Let us now consider QCD in the finite volume V where \AB, W) and \CD,W) 
are no longer eigenstates of the Hamiltonian. True eigenvalues are shifted from W 
to Wi = W + OiV^ 1 ) (i = 1,2). By the diagonalization method in lattice QCD, 
it is relatively easy to determine W\ and Wi- With these values Liischer's finite 
volume formula gives two conditions, which, however, are insufficient to determine 
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three observables, two scattering phase shifts 8\, 5j and one mixing angle 9. We here 
explain a new approach proposed in Refs.60),61) to overcome this difficulty. Let us 
consider the NBS wave functions at two different values of energy, W\ and W2, in 
the finite volume: 

Vx (x,pf)e- Wit = (0\T{ l f Xl (x + r,t) VX2 (r,t)}\W i ), i = 1,2. (7-7) 

where X{= X 1 X 2 ) = AB or CD with p AB = k or p CD = q. We then define the 
coupled channel non-local potentials from the coupled channel Schrodinger equation 
as 

<Px(x,p?) = Y, d3 V U x,v(x,y) <p Y {y,pJ) (7-8) 

Y J 

for i = 1,2 where the reduced mass is defined by 1/nx = V m ^i + V m A" 2 - I n the 
leading order of the velocity expansion, we have 

<Px(x,p?) = Y / V x ,y(x) <f Y (x,pJ) (7-9) 

Y 

These equations for i = 1, 2 can be solved as 

( V A b,ab{x) V A b,cd(x) \ = ( K AB (x,k 1 ) 
V Vcd,ab{x) V C d,cd(x) J \ KcDfaQl) 

x / <^4s(a;,fci) 
\ ifCDix^q-i) 

Once we obtain the coupled channel local potentials Vx,y( x ), we solve the cou- 
pled channel Schrodinger equation in infinite volume with some appropriate bound- 
ary condition such that the incoming wave has a definite £ and consists of the AB 
state only, in order to extract three observables for each £ (5j(W), Sj(W) and 9g(W)) 
at all values of W. Of course, since Vx,y is the leading order approximation in the ve- 
locity expansion of Ux,y( x , y), results for three observables 5\{W), 5j (W) and Q(W) 
at W 7^ Wi, W2 are also approximate ones and might be different from the exact 
values. By performing an additional extraction of V x ,y( x ) at (W3, W4) / (W±, W2), 
we can test how good the leading order approximation is. 

The method considered above can be generalized to inelastic scattering where a 
number of particles is not conserved such that A+B — > A+B and A+B — > A+B+C. 
See Ref. 61) for more details. 

7.2. Coupled- channel potentials in (S, I) = (—2, 0) channel 

As an application of the method in the previous subsection, let us consider BB 
potentials for the S = —2 and 1 = system, i.e., the coupled AA-NESU system. 
It is interesting to investigate this system since it involves the flavor single state, 
which is free from the Pauli blocking of quark degrees of freedom at short distance. 
Since mass differences of these BB systems are quite small, all 3 states appear in 



(p 



A Y2 



2»x 



He, 



K x (x,p?) = 



X\2 



K AB (x,k 2 ) 
K CD (x,q 2 ) 

ip AB (x,k 2 ) 
ip C D{x,q 2 ) 



(7-10) 
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Fig. 19. The coupled channel potential matrix from the NBS wave function for Set 1. The vertical 
axis is the potential strength in units of MeV, while the horizontal axis is the relative distance 
between two baryons in units of fm. 



NBS wave functions. The coupled channel method in the previous subsection can 
be applied to treat such complicated systems. 

At the leading order of the velocity expansion, the coupled channel 3x3 potential 
matrix in this case is given by 

3 r i-l 

i=i 

where i is a label for energy W% and X,Y = AA, NS or UU. Here the last factor is 
the inverse of the 3x3 matrix ip^ir, k A ) with indices i and A. 

Gauge configurations generated by CP-PACS/JLQCD Collaborations on a 16 3 x 
32 lattice at a ~ 0.12 fm ( therefore L ~ 1.9 fm) in 2+1-flavor full QCD simulations 
are employed to calculate the coupled channel potentials at three different values of 
the light quark mass with the fixed bare strange quark mass. 62 ) Quark propagators 
are calculated with the spatial wall source at to with the Dirichlet boundary condition 
in time at t = to + 16. Corresponding hadron masses are given in Table IV. 

The coupled channel potential matrix Vab,cd from the NBS wave function for 
Set 1 is shown in Fig. 19. All diagonal components of the potential matrix have a 
repulsion at short distance. The strength of the repulsion in each channel, however, 
varies, reflecting properties of its main component in the irreducible representation 
of the flavor SU(3). In particular, the SJS potential has the strongest repulsive core 
of these three components. It is important to note that off-diagonal parts of the 
potential matrix satisfy the hermiticity relation Vab,GD = Vcd,ab within statistical 
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Set 1 
Set 2 
Set 3 


875(1) 
749(1) 
661(1) 


916(1) 
828(1) 
768(1) 


1806(3) 
1616(3) 
1482(3) 


1835(3) 
1671(2) 
1557(3) 


1841(3) 
1685(2) 
1576(3) 


1867(2) 
1734(2) 
1640(3) 



Table IV. Hadron masses in units of MeV and number of configurations for each set adopted in 
Sec.7.2. 



Lattice QCD approach to Nuclear Physics 



33 




r[fm] 



25 0.5 0.75 

<m 



# + + 



s ° - 



0.25 0.5 0.75 

r[lm] 



0.25 0.5 0.75 



r[fm] 



^ Set! - 

11 SU(3) - 

'» V _ 



m 



0.25 0.5 0.75 1 1.25 1.5 



0.25 0.5 0.75 1 1.25 1.5 

r[fm] 



0.2 0.4 0.6 0.8 



Fig. 20. Transition potentials in the flavor SU(3) IR basis. Red, blue and green symbols correspond 
to results of Setl, Set2 and Set3, respectively. The result of the flavor SU(3) symmetric limit 
at the same strange quark mass is also plotted with brown symbols 52). 



errors. In addition the off-diagonal parts of AA to NS transition Vaa,ne is much 
smaller than the other two off-diagonal potentials, Vaa,ee and Vns,ee- 

In order to compare the results of the potential matrix calculated in three con- 
figuration sets, the potentials from the particle basis are transformed by the unitary 
rotation U to those in the flavor SU(3) irreducible representation (IR) basis as 

/ ^1,1 ^1,8 V!,27 \ 

V IR = U^VU = V 8 ,i Vs, 8 V 8)27 ■ (7-12) 

\ Va7,l ^27,8 ^27,27 / 

The potential matrix in the IR basis give a good measure of flavor SU(3) breaking 
effects since it is diagonal in the SU(3) symmetric limit. 

In Fig. 20, the results of the potential matrix in the IR basis are compared among 
different configuration sets, together with the one in the flavor SU(3) symmetric 
limit. As the pion mass decreases, the attraction in Via potential increases in whole 
range. While the Vs t s potential in the flavor SU(3) limit deviates from others, we 
do not observe a clear pion mass dependence of the potentials among these flavor 
breaking cases. In the ^7,27 potential, we observe the growth of attraction range 
of potential and the enhancement of repulsive core. The V\ t 27 and Vs,27 transition 
potentials are consistent with zero within statistical errors. On the other hand, it is 
noteworthy that the flavor SU(3) symmetry breaking effect becomes manifest in the 
Vi 8 transition potential. 

§8. Three-Nucleon Forces 

In this section, we expand the scope of our studies from two-nucleon (2N) sys- 
tems to vl-body nucleon systems. Generally speaking, there could exist not only 
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2-body forces but also 3, 4, • • • ^4-body forces in such systems. In particular, the de- 
termination of three-nucleon forces (3NF) attracts a great deal of interest, since it has 
been revealed that 3NF play an important and nontrivial role in various phenomena. 
Some examples include (a) binding energies of light nuclei, 63 ) (b) deuteron-proton 
elastic scattering experiments, 64 ) (c) the anomaly in the oxygen isotopes near the 
neutron drip-line, 65 ) and (d) the nuclear equation of state (EoS) at high density 
relevant to the physics of neutron stars. 66 ) Universal short-range repulsion for three 
baryons (nucleons and hyperons) has also been suggested in relation to the maximum 
mass of neutron stars with hyperon core. 67 )' 68 ) 

Despite of its phenomenological importance, microscopic understanding of 3NF 
is still limited. Pioneered by Fujita and Miyazawa, 69 ) the long range part of 3NF has 
been modeled by the two-pion exchange (27rE), 70 ) particularly with the Z\-resonance 
excitation. This 2-7rE-3NF component is known to have an attractive nature at 
long distance. An additional repulsive component of 3NF at short distance is often 
introduced in a purely phenomenological way. 71 ) An approach based on the chiral 
effective field theory (EFT) is quite useful to classify the two-, three- and more- 
nucleon forces and has been studied intensively. 72 ) . A completely different approach 
based on holographic QCD is recently proposed, which obtains repulsive 3NF at short 
distance. 73 ) 

To go beyond phenomenology, it is most desirable to determine 3NF directly 
from the fundamental degrees of freedom (DoF), the quarks and the gluons, on the 
basis of quantum chromodynamics (QCD). In this section, we present an exploratory 
study of first-principle lattice QCD calculation of 3NF in the quantum numbers of 
(I, J p ) = (1/2, 1/2+) (the triton channel). For details of this study, see Ref. 74). 

8.1. Formalism 

We consider the NBS wave function ip 3 ^(r,p) extracted from the six-point cor- 
relator as 

G 3N (r,p,t-t ) = -L J2(0\(N( Xl )N(x 2 )N(x 3 ))(t) (N'N'N')(t )\0), (8-1) 

R 



A m ip 3N (r,p)e- E3N(t - t0 \A m = (E 3N \(N'N'N'){0)\0), 



t>to 

^3Jv(r,p) = (0\N(x 1 )N(x 2 )N(x 3 )(0)\E 3N ), (8-2) 

where E 3 n and \E 3 n) denote the energy and the state vector of the 3N ground 
state, respectively, N (N 1 ) the nucleon operator in the sink (source), and R = 
(xi + X2 + x 3 )/3, r = X\ — X2, p = x 3 — {x\ + X2)/2 the Jacobi coordinates. Wc 
consider the following Schrodinger equation of the 3N system with the derivative 
expansion of the potentials, 



" 27T V ' " 27T V ' + £ V ™{Tij) + V 3 NF(r, p) 
Mr H-p i< - 



*p3N(r,p) = E 3N ip 3N (r,p), 

(8-3) 

where V^N^y) with = x, L — xj denotes the 2NF between (i,j)-pair, V 3 NF{r,p) 
the 3NF, \i T = mjv/2, \i p = 2m^/2> the reduced masses. 3NF can be determined as 
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follows. We first calculate "03Jv( r >P) and obtain the total potential of the 3N system 
through Eq. (8-3). Once we obtain all necessary V^/v^i?) by performing (separate) 
lattice simulations for genuine 2N systems, we can extract Vsnf(t, p) by subtracting 
Si<j ^2N( r ij) from the total potential. The extension to four- and more-nucleon 
forces can be immediately understood. Note that potentials determined in this way 
reproduce the energy of the system by construction. 

An important remark is that 3NF are always determined in combination with 
2NF, and 3NF alone do not make too much sense. Therefore a comparison between 
lattice 3NF and phenomenological 3NF can be done only at a qualitative level. 
Rather, our purpose is to determine two-, three-, (more-) nucleon forces systemati- 
cally, and to provide them consistent set. 

One of the difficulties in the 3NF study from lattice QCD is that computational 
costs become exceptionally enormous. Since there are 9 valence quarks, the DoF of 
color and spinor are significantly enlarged. In addition, the number of diagrams in 
the Wick contraction tends to diverge with a factor of N u \ x N^\, where N u (Nj) are 
numbers of up (down) quarks in the system. We here develop several techniques to 
reduce these computational costs. We first take an advantage of symmetries (such as 
isospin symmetry) to reduce the number of Wick contractions. Second, we utilize a 
freedom for the choice of a nucleon interpolating operator. In particular, a potential 
is independent of the choice of a nucleon operator at the source, TV 7 , and we employ 
the non-relativistic operator as N' = N nr = e a bc{Qa C^^P nr qb)P n rqc with P nr = 
(1 + 74)/2, which reduces the spinor DoF. Similar techniques are (independently) 
developed in Ref. 75). On the other hand, a potential depends on the choice of a 
nucleon operator at the sink, N. As discussed in Sec. 2.4, choosing N corresponds 
to choosing the "scheme" to calculate nuclear forces. Note that physical observables 
calculated from these different potentials such as phase shifts and binding energies 
are unique. In order to determine 3NF and 2NF in the same "scheme", we employ 
the same sink operator N = N stc [ = e a b c (q^ 'C75 in the 3NF study, as employed 
in 2NF. Recall that the choice of N = N st d is shown to have good convergence in 
the derivative expansion in Sec. 4.3, and can be considered to be a good "scheme" 
for lattice nuclear forces. 

We next consider the geometry of the 3N. Since the spacial DoF for general 
(r, p) is too large, it is necessary to find an efficient way to restrict the geometry. In 
this exploratory study, we propose to use the "linear setup" with p = 0, with which 
3N are aligned linearly with equal spacings of r2 = \r\/2. In this setup, the third 
nucleon is attached to (1, 2)-nucleon pair with only S-wave. Considering the total 
3N quantum numbers of (I, J p ) = (1/2, l/2 + ), the triton channel, the wave function 
can be completely spanned by only three bases, which can be labeled by the quantum 
numbers of (l,2)-pair as 1 So, 3 Si, 3 Di. Therefore, the Schrodinger equation leads 
to the 3x3 coupled channel equations with the bases of ^>ig , ^35^, ipa Dl . The 
reduction of the dimension of bases is expected to improve the S/N as well. It is 
worth mentioning that considering the linear setup is not an approximation: Among 
various geometric components of the wave function of the ground state, we calculate 
the (exact) linear setup component as a convenient choice to study 3NF. While we 
can access only a part of 3NF from it, we plan to extend the calculation to more 
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Fig. 21. 3N wave functions at (t — to) /a = 8. 
Circle (red), triangle (blue), square (green) 
points denote ips, ipM, V^Dii respectively. 
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Fig. 22. The effective scalar-isoscalar 3NF in 
the triton channel with the linear setup ob- 
tained at (t — to) /a — 8. r2 is the distance 
between the center and edge in the linear 
setup. 



general geometries step by step, toward the complete determination of the full 3NF. 

Finally, we emphasize that 3NF study requires the precise determination of 2NF. 
This is not surprising because the interactions in 3N systems are mostly dominated by 
2NF, and thus small uncertainties in 2NF could easily obscure the signal of 3NF. Note 
that we generally need precise 2NF in both of parity-even and parity-odd channels, 
since a 2N-pair inside the 3N system could be either of positive or negative parity. 
On this point, we have shown that 2NF in parity-even channel can be determined 
with good precision in lattice QCD. On the other hand, the determination of 2NF 
in parity-odd channel is much more difficult. While the formulation for parity-odd 
2NF is developed in Sec. 4.5, the results are found to suffer from larger statistical 
errors than parity-even 2NF. This is considered to be a general tendency, since one 
has to inject a non-zero momentum in the parity-odd 2NF study. Therefore, in the 
3NF study, it is essential to suppress the uncertainties originated from parity-odd 
2NF. 

In order to address this issue, we propose to consider the following channel, 74 -* 



•4) 



which is anti-symmetric in spin/isospin spaces for any 2N-pair. Combined with 
the Pauli-principle, it is automatically guaranteed that any 2N-pair couples with 
even parity only. Therefore, parity-odd 2NF vanish in (ips\H\ip3N)j where H is the 
Hamiltonian of the 3N system, and we can extract 3NF unambiguously without 
referring to parity-odd 2NF. Note that no assumption on the choice of 3D-geometry 
of r, p is imposed in this argument, and we can take an advantage of this feature 
for future 3NF calculations with various 3N geometries. 

8.2. Numerical results 

We employ Nf = 2 dynamical configurations with mean field improved clover 
fermion and RG-improved gauge action generated by CP-PACS Collaboration. 39 ) We 
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use 598 configurations at f3 = 1.95 and the lattice spacing of a -1 = 1.269(14) GeV, 
and the lattice size of V = L 3 x T = 16 3 x 32 corresponds to (2.5 fm) 3 box in physical 
spacial size. For u, d quark masses, we take the hopping parameter at the unitary 
point as k U( i = 0.13750, which corresponds to m n = 1.13 GeV, = 2.15 GeV 
and niA = 2.31 GeV. We use the wall quark source with Coulomb gauge fixing, and 
periodic (Dirichlet) boundary condition is imposed in spacial (temporal) direction. In 
order to enhance the statistics, we perform the measurement at 32 source time slices 
for each configuration, and the forward and backward propagations are averaged. 
The results from both of total angular momentum J z = ±1/2 are averaged as well. 
Due to the enormous computational cost, we can perform the simulations only at 
a few sink time slices. Looking for the range of sink time where the ground state 
saturation is achieved, we carry out preparatory simulations for effective 2NF in the 
3N system 74 ) in the triton channel at 2 < (t — to) /a < 11, and find that the results 
are consistent with each other as long as (t — to) /a > 7. 74 ) Being on the safer side, 
we perform linear setup calculations of 3NF at (t — to) /a = 8 and 9. We perform 
the simulation at eleven values of the spacial distance r^. 

In Fig. 21, we plot the radial part of each wave function of ips = (—^s + 
V'35 1 )/\/2, i>M = (V^So + V' 3 5i)/v / 2 and obtained at (t — to) /a = 8. Here, we 
normalize the wave functions by the central value of ips( r 2 = 0). What is noteworthy 
is that the wave functions are obtained with good precision, which is quite nontrivial 
for the 3N system. We observe that i/js overwhelms the wave function, and tpM, 
ipso! are much smaller by one to two orders of magnitude. This indicates that 
higher partial wave components in tps are also strongly suppressed, and the wave 
function is completely dominated by the component with which all three nucleons 
are in S-wave in this lattice setup. 

We determine 3NF by subtracting 2NF from total potentials in the 3N sys- 
tem. As discussed before, we have only one channel, {ips\H\ip3N) = {4's\H\'4's) + 
(ifjs\H\tpM) + (V'si-HlV ; 3£) 1 ), which is free from parity-odd 2NF. Correspondingly, 
we can determine one type of spin/isospin functional form for 3NF. In this study, 
3NF are effectively represented in a scalar-isoscalar functional form. This form 
is an efficient representation, since overwhelms the wave function and thus 
\{ips\V3NF\ips)\ > \{iPs\V3Nf\iI>m)\,\(' 1 I ; s\V3nf\iIJ3d 1 )\ & expected. Note also that 
the scalar-isoscalar functional form is often employed for the short-range part of 
3NF in phenomenological models. 71 ) 

In Fig. 22, we plot the results for the effective scalar-isoscalar 3NF at (t — to)/a = 
8. We here include ^-independent shift by energy, Se — 5 MeV, which is determined 
by long-range behavior of potentials . 74 ) While Se suffers from <J 10 MeV systematic 
error, it does not affect the following discussions much, since Se merely serves as an 
overall offset. In order to check the dependence on the sink time slice, we compare 
3NF from (t — to) /a = 8 and 9. While the results with (t — to) /a = 9 suffer from 
quite large errors, they are consistent with each other within statistical fluctuations. 

Fig. 22 shows that 3NF are small at the long distance region of r2- This is in 
accordance with the suppression of 27rE-3NF for the heavy pion. At short distance, 
however, an indication of repulsive 3NF is observed. Note that a repulsive short- 
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range 3NF is phenomenologically required to explain the properties of high density 
matter. Since multi-meson exchanges are strongly suppressed for the large quark 
mass, the origin of this short-range 3NF may be attributed to the quark and gluon 
dynamics directly. In fact, we recall that the short-range repulsive (or attractive) 
cores in the generalized two-baryon potentials in the flavor SU(3) limit discussed in 
Sec. 6 are consistent with the Pauli exclusion principle in the quark level. 45 )' 52 ) In 
this context, it is intuitive to expect that the 3N system is subject to extra Pauli 
repulsion effect, which could be an origin of the observed short-range repulsive 3NF. 
Further investigation along this line is certainly an interesting subject in future. 

We remark here that the quark mass dependence of 3NF is certainly an impor- 
tant issue, since the lattice simulations are carried out only at single large quark 
mass. In the case of 2NF, short-range cores have the enhanced strength and broaden 
range by decreasing the quark mass. 18 ) We therefore would expect a significant 
quark mass dependence exist in short-range 3NF as well. In addition, long-range 
2-7rE-3NF will emerge at lighter quark masses, in particular, at the physical point. 
Quantitative investigation through lattice simulations with lighter quark masses are 
currently underway. 

§9. Meson-baryon interactions 

Since our potential method can be naturally extended to meson-baryon systems 
as well, we consider meson-baryon interactions in this section. The first application 
is a study on kaon-nucleon (KN) interactions in the I(J P ) = 0(l/2~) and l(l/2~) 
channels. The elastic KN scattering allows us to study the origin of "non-resonant" 
nuclear forces, since kaon contains us quarks, and these quarks do not annihilate 
in non-strange nucleons. Therefore, the direct productions of conventional baryon 
resonances are ruled out. Also, the KN systems in the I(J P ) = 0(l/2~) and l(l/2~) 
channels may be relevant for a possible exotic state 6> + , 76 ) whose existence is still 
controversial. 

It is important to emphasize that the one-pion exchange is absent in the KN 
systems, so that the short- and mid-range interactions dominate the elastic KN 
scattering. Theoretical studies of the KN interactions so far have been carried 
out by constituent quark models and meson-exchange models. In both models, it 
was found that genuine quark-gluon dynamics become important to describe the 
empirical scattering phase shifts. 77 ) 

To investigate the KN potentials in 2+1 flavor full QCD, we have utilized the 
gauge configurations of JLDG(Japan Lattice Data Grid) /ILDG (International Lattice 
Data Grid) generated by PACS-CS Collaboration on a 32 3 x 64 lattice. 78 ) The 
renormalization group improved Iwasaki gauge action and non-perturbatively 0(a) 
improved Wilson quark action are used at j3 = 1.90, which corresponds to the lattice 
spacing a = 0.09 fm determined from 7r, K and fl masses. The physical size of 
the lattice is about (2.9 fm) 3 and the the hopping parameters are taken to be k u = 
Kd = 0.1370 and k s = 0.1364. In the present simulation, we adopt the wall source 
located at to with the Dirichlet boundary condition at time slice t = to + 32 in the 
temporal direction and the periodic boundary condition in each spatial direction. 
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Fig. 23. The LO S-wave central potentials for the KN states in the 1 = (left) and the I 
channels (right). 



The Coulomb gauge fixing is employed at t = to. The number of gauge configurations 
used in the simulation is 399. With this setup, we obtain m n = 705(2), m# = 793(2) 
and M N = 1590(8) MeV. 79 ) 

Fig. 23 shows the LO S-wave potential V(r) for the KN state in the 1 = (left) 
and 1 = 1 (right). The repulsive interactions are observed at short distance in both 
channels, while the attractive well appears at the mid-distance (0.4 < r < 0.8 fm) in 
the 1 = channel, which is not found in the constituent quark model of hadrons. 80 ) 
These results indicate that there are no bound states in I(J n ) = 0(l/2~) and l(l/2~) 
states at ~ 705 MeV. 

The strong repulsions near origin in the 1 = 1 channel can be expected by the 
quark Pauli blocking effects. This was first pointed out by Machida and Namiki 81 ) 
for the meson-baryon systems. In the 1 = 1 KN (K + p) state whose configuration is 
K + p ~ (us) (uud), one of the u-quarks cannot be in the S-state. The strong repulsion 
at short distance in this channel found in our simulations suggests a manifestation of 
the quark Pauli blocking. In addition, the repulsive interactions in the S-wave 1 = 1 
KN state can be expected much stronger than that of the 1 = KN state in the 
constituent quark model. Our simulation shows that the repulsion at short distance 
for the KN potential becomes significantly smaller in the 1 = channel than 1 = 1 
channel. This again confirms the expectation from the quark Pauli blocking effects. 

Kawanai and Sasaki investigate charmonium-nucleon (cc-N) interactions using 
the potential method. Since charmonia do not share the same quark flavor with the 
nucleon, the cc-N interaction is mainly induced by the genuine QCD effect of multi- 
gluon exchanges and does not manifest the repulsive core near origin. In Ref. 82), 
the charmonium-nucleon potentials are calculated in quenched QCD on 16 3 x 48 and 
32 3 x 48 lattices at a ~ 0.94 fm at three different values of the light quark mass 
corresponding to (m^,m N ) ~ (640, 1430), (720, 1520), (870, 1700) in unit of MeV 
and one fixed value of the charm quark mass corresponding to m„ c ~ 2920 MeV 
and mj/q, ~ 3000 MeV. They have found that the effective central cc-N potentials 
clearly exhibit entirely attractive interactions without any repulsion at all distances. 
Absence of the short range repulsion (the repulsive core) is related to absence of the 
Pauli exclusion between the heavy quarkonium and the light hadron. 
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§10. Conclusion 

In this report, we reviewed the basic notion of the HAL QCD method for the 
baryon-baryon (BB) potential and its field-theoretical derivation from the equal-time 
Nambu-Bethe-Salpeter (NBS) wave function in QCD. The potential U(x,y) (or the 
integral kernel of the Schrodinger type equation) in the HAL QCD method has three 
characteristic features: (i) non-local in relative coordinate, (ii) energy independent, 
and (hi) scheme dependent. Each of these features has been discussed in detail in 
this report. By construction, this potential correctly reproduces the scattering phase 
shift defined from the S-matrix in QCD below the inelastic threshold. 

One can construct U(x,y), once all the NBS wave functions for scattering en- 
ergies below the inelastic threshold are obtained. In lattice QCD simulations in a 
finite box, however, it is more practical to adopt the velocity (derivative) expansion 
of U(x, y) by its non-locality and determine the local potentials V(x) order by order. 
This is also in conformity with phenomenological potentials widely used in nuclear 
physics: An advantage of the HAL QCD method is that one can check the accuracy 
of this velocity expansion by changing the scattering energies on the lattice. To avoid 
the well-known problem of exponential error-growth in the temporal correlation of 
multi hadrons, we have introduced a time-dependent HAL QCD method on the basis 
of the (imaginary) time Schrodinger type equation. Due to this improved method, 
we could achieve better construction of the potential as demonstrated in this report. 

The leading order (LO) terms of the velocity expansion correspond to the cen- 
tral and tensor potentials: Those for the nucleon-nucleon (NN), hyperon-nucleon 
(YN) and hyperon-hyperon (YY) interactions have been investigated in full QCD 
simulations, some of which are recapitulated in this report. The next-to- leading or- 
der (NLO) term is the spin-orbit potential: By introducing finite momentum to the 
nucleons, we could extract the NN spin-orbit force for the first time. The origin of 
the repulsive core of the NN interaction has been also investigated by extending the 
SU(2)-flavor to degenerate SU(3)-flavor. The role of the Pauli principle in the quark 
level to describe the short range part of the interaction becomes clear. In particular, 
there arises a short range "attractive" core in the flavor singlet channel; we found 
that it is strong enough to form a bound state, i^-dibaryon, in the SU(3) limit. 

The HAL QCD method can be extended to the case beyond the inelastic thresh- 
old. This is necessary to treat the YN and YY interactions with SU(3)-flavor sym- 
metry breaking. We have presented its application to (S, I) = (—2, 0) system and 
derived the coupled channel potentials among AA, NS and EE. The HAL QCD 
method is also applied to the three-nucleon force relevant for the extra binding of fi- 
nite nuclei and also for the maximum mass of neutron star, and to the meson-baryon 
interactions relevant for the meson-baryon resonances and the pentaquark. 

So far, our full QCD simulations of the BB interactions are performed at non- 
zero lattice spacing on a finite volume with relatively large quark masses. We there- 
fore need careful studies of systematic errors on finite volume effect, quark mass 
dependence and the lattice spacing effect. Among others, the most important di- 
rection is to carry out (2+l)-flavor simulations on a large volume (e.g. L = 6 — 9 
fm) at physical quark mass (m, = 135 MeV) to extract the realistic BB and BBB 
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potentials. Such simulations are planned at 10 PFlops "K computer" in Advanced 
Institute for Computational Science (AICS), RIKEN. 

If it turns out that the program described in this paper indeed works in lattice 
QCD at the physical quark mass, it would be a major step toward the understanding 
of atomic nuclei and neutron stars from the fundamental law of the strong interaction, 
the quantum chromodynamics. 
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